Programming Language:
    Complex Numbers | Linear Algebra | Fit Algorithms | Interpolation | Root Finding | ODE | FFT | Special Functions | Integration

    LINEAR ALGEBRA

    - standard linear algebra algorithms for double and complex values in Xi-

    Define an example vector and matrix and transpose the matrix

    (  1)>A={ {1,3}, {2,4} }; B=transpose(A);
    (  2)>b={1,5};
    
    To perform a matrix multiplication use the "#" operator
    (  3)>print(A#A);
    <intarr>
     7 15 
    10 22
    
    The linear equation A*x=b will be solved by
    (  4)>print(linear_solve(A,b));
    <dblarr>
     5.5 -1.5
    
    The command
    (  5)>print(invert(A));
    <dblarr>
      -2  1.5 
       1 -0.5 
    
    gives the inverse of the matrix A. Furthermore to compute the determinant and to estimate the reciprocal condition numbers (first to the 1-norm, second to the infinity-norm) type
    (  6)>print(determ(A),rcond(A),rcond_inf(A));
    <double> -2
    <double> 0.047619
    <double> 0.047619
    
    Just as a reminder: rcond(A)=1 / ( norm(A) * norm(inv(A)) ). The eigenvalues and eigenvector can be obtained by
    (  7)>print(eigenvalue(A),eigenvector(A));
    <cpxarr>
    (  5.3722813, 0) (-0.37228132, 0)
    <dblarr>
     0.56576746 -0.90937671 
     0.82456484  0.41597356
    
    The eigenvectors are stored in the columns of the result matrix. The command
    (  8)>[LU, pivot]=lu(A);
    (  9)>print(LU, pivot);
    <dblarr>
      2   4 
    0.5   1 
    <intarr>
    2 2
    
    performs a LU decomposition of the matrix A. This is the first example where a function has more than one (here two) return arguments. The first return value LU contains the factors L and U from the factorization A = P*L*U; the unit diagonal elements of L are not stored. The second one represents the pivot indices; row i of the matrix A was interchanged with row pivot[i]. Proceed this computation you can evalute the inverse of A by
    ( 10)>inv=lu_invert(LU, pivot);
    
    If you do not need the pivot indices, call the function by
    ( 11)>LU=lu(A);
    
    A QR decomposition can be computed by,
    ( 12)>[t, R]=qr(A);
    <dblarr>
    1.4472136         0
    <dblarr>
      -2.236068  -4.9193496 
     0.61803399 -0.89442719
    
    where the first return array has the scalar factors of the elementary reflectors as entries and the second value contains the Q and R matrix (for further details see the reference guide). Another supported decomposition is the SVD.
    ( 13)>[w,u,vt]=svd(A);
    
    Here w is a two element vector of singular values, u a (2 by 2) matrix and vt a (2 by 2) orthogonal matrix. A can be rewritten as
    ( 14)>W={ {w[0],0} , {0,w[1]} };
    ( 15)>print(u#W#vt);               /* u#W#vt equal to A */
    <dblarr>
    1 3 
    2 4
    
    Uses "back substitution" to solve a linear equation A*x=b
    ( 16)>x=svd_back_sub(w,u,vt,b);
    
    A typical use of the functions svd and svd_back_sub is to zero the very small w[i]'s before executing the back substitution to make the result stable against small changes in b. (for further details have a look in a standard book on matrix numerics).


    Rechts Index Index Index Linls © 1995 by Bodo Junglas, Klaus Spanderen and Fabian Weis
    - Last revised: April 23 1996